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Rock physics and geophysics for unconventional resource, multi-component 
seismic, quantitative interpretation 
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An extension of a previously developed, rock physics, model is made that quantifies the relationship between 
the ductile fraction of a brittle/ductile binary mixture and the isotropic seismic reflection response. Making 
a weak scattering (Born) approximation and plane wave (eikonal) approximation, with a subsequent ordering 
according to the smallness of the angle of incidence, a linear singular value decomposition analysis is done to 
understand the stack weightings, number of stacks, and the type of stacks that will optimally estimate the 
two fundamental rock physics parameters. It is concluded that the full PP stack and the "full" PS stack are 
the two optimal stacks needed to estimate the two rock physics parameters. They dominate over both the 
second order AVO "gradient" stack and the higher order (4th order) PP stack. 



I. INTRODUCTION 

The developing commercial significance of unconven- 
tional shale reservoirs is leading to the need to be able to 
remotely determine the ability to effectively fracture the 
reservoir. This paper will establish the theory and practi- 
cality of optimally estimating the ductile fraction from an 
isotropic analysis of surface conventional and converted 
wave seismic. This property of a binary ductile/brittle 
mixture has been shown to be the key property in de- 
termining the geomechanical fracturing response of an 
unconventional reservoir^. This is, most likely, because 
of the balance between the "bumpy road" friction of the 
fracture, due to the structurally competent brittle mem- 
ber, and the viscous friction, due to the ductile member. 
This is not the subject of this paper, but is the topic 
of our ongoing research into the statistical mechanics of 
fracture joint friction. 

Because of unrelated physics, the same property, duc- 
tile fraction, is one of two important order parameters 
for the linear, isotropic, elastic response of binary mix- 
tures of a structurally competent member (high coordina- 
tion number) and a structurally less competent member 
(lower coordination number). A very important impli- 
cation of this bicritical model is that the state is only 
two dimensional. The expectation, and practical reality 
(as demonstrated by analysis of well log data) is that 
the isotropic properties will reduce to a surface in the 
three dimensional density, compressional velocity, shear 
velocity (i.e., p,Vp,Vs) space. Furthermore, this surface 
will be orthogonal to the Vp-v^ plane. This remark- 
able property is captured by the floating grain modeP^l 
which has two state variables given by the floating grain 
fraction, ff, and the compaction state as specified by 
1 — exp(— Pe/Po)i where P^ is the effective stress and 
Pq is a reference value of effective stress. Two phase 
transitions points at critical values in the radius ratio 
{RRc = 4) and the fraction of small grains {VFc — 0.45) 
were demonstrated, as well as two critical scalings of the 
porosity about a critical point of about 429P1. 

This theory was developed for a binary mixture of 
brittle spheres of two different sizes. Recognizing that 



the large spheres are the structurally competent member 
and the small spheres are the structurally less compe- 
tent member, we generalize this theory in Sec. [TTj The 
floating grain fraction is replaced by a general geome- 
try parameter, ^, which in the case of shales is shown 
to be proportional to the \/Yd, where fii is the ductile 
fraction. The geometry parameter captures the fabric of 
the mixture such as the sorting or ductile fraction, while 
the composition parameter captures the compaction, di- 
agenesis, and/or mineral substitution of the mixture. An 
important additional implication of the bicritical model 
is a fundamental self similarity and the associated scaling 
relationships of physical quantities such as coordination 
numbers, capture fractions, and elastic moduli. It also 
implies the same critical scaling for both Vp and v^ be- 
cause they have the same units. Therefore the surface in 
{p,Vp,Vs) space must be orthogonal to the Vp-Vg plane. 

We emphasize the serendipity of the fact that the duc- 
tile fraction is the coordinate of influence of both the 
linear elastic response (geophysical) and the nonlinear 
inelastic response (geomechanical). For the former, the 
ductile material is adding density without much struc- 
tural rigidity, that is elastic moduli. For the latter, it is 
increasing the importance of the viscous joint friction. 

Given this rock physics model, this paper examines 
its implication on the geophysical detectability of ductile 
fraction in Sec. 



Ill 



Several questions have been the sub- 
ject of much debate within the geophysical community. 
How many stacks should be used in "prestack" analysis? 
What should those stacks be? What is the relative value 
of AVO versus converted wave analysis? What are the 
quantities that should be inverted for, relative (reflectiv- 
ity) versus absolute (impedances)? Finally, what are the 
"attributes" that best predict reservoir performance? 

We present a straight forward analytic theory in Sec. 



IV and subsequent analysis that answers all of these ques- 



tions. It is a linear singular value decomposition analy- 
sis of the relationship between the two fundamental rock 
physics parameters (C and ^) and the seismic reflectiv- 
ities (PP and PS) as functions of angle of incidence, 6. 
This analysis is done assuming a weak scattering (Born) 
approximation and plane wave assumption (eikonal). It 



also orders the SVD using the angle 9. The conclusion is 
that the PP full stack and the PS "full" (linear weighted 
with or offset) stacks are optimal in the estimation of C, 
and ^, respectively. They are of zeroth and first order in 
0, respectively. Conventional AVO "gradient" stacks and 
large angle PP response conventionally used to estimate 
density are of higher order in (second and fourth order, 
respectively). It should be noted that these stacks are 
average reflectivities or relative quantities. Linear com- 
binations of these two important stacks (normally just 
the full PP stack for C and the "full" PS stack for are 
the best "attributes" . 



II. ROCK PHYSICS 

We first recognize that we are dealing with a binary 
mixture of a ductile and a brittle member, where the 
latter is more structurally competent than the former. 
We now take inspiration from the floating grain modeH 
This model is based on two fundamental parameters - the 
floating grain fraction parameterized by ^ = If I ffc and 
the compaction parameterized by C = 1 — exp(— Pe/-Fb)i 
where ff is the floating grain fraction, ff^ is the maxi- 
mum or critical floating grain fraction, Pg is the effective 
stress, and Pq is a reference effective stress. The model 
respects fluid substitution and leads to local linear cor- 
relations of the form 



Vp — -^vp 



Byp C, + Cvp £^ zLa. 



' vpj 



^ A^ + B^Vp + C^^ ± I 



and 



(1) 

(2) 
(3) 



The second relationship can be rewritten two ways, given 
Ps and pf and the definition p = (j) pf + {1 — (f>)ps. 



P = Ap 



■c- 



BpVp 



— ^, and 

f Cp^ ±<7p. 



(4) 
(5) 



The first identifies the two critical exponents, n^ and 
n^, and the critical porosity, 4>c- The second is just a 
convenient expression to compare to log data of shales. 

When we examine shales from many different wells and 
plays, we get the results shown in Fig. [T] It is important 
to note the strong linear correlation in the Vp-Vg plane of 
Fig. Ilk, and the systematic shift in the p-Vp correlation 
with the ductile fraction, /d, in Fig. [l]3. Inspired by 
the floating grain model, we generalize C to {fd/ fdc)"''', 
where fd is the ductile fraction, fdc is the maximum or 
critical ductile fraction, and Ud is the critical exponent 
for the ductile fraction. Given the range of the data, we 
assume the line in Fig. [Tb shows the variation in the 
p-Vp as f goes from to f for C = I, and we assume that 
the minimum value of Vp is 9500 ft/s when ^ = ^ = 0. 
A regression to this extended rock physics model leads 
to Ayp = 9500 ft/s, Byp = 8500 ft/s, C^p = -4500 ft/s. 
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FIG. 1. Well log data supporting rock physics model. Points 
are blocked well log data colored according to the ductile 
fraction, fd- Also shown are the directions of increasing ( 
(constant ^) as the red arrow, and increasing ^ (constant () 
as the green arrow. Values are normalized according to the 
equation x = {x — x,nin)/ixmax — Xmin), where minvp = 8000 
ft/s, max Up = 18000 ft/s, rninu^ = 3800 ft/s, maxu^ = 11000 
ft/s, minp — 2.1 gm/cc, maxp — 2.8 gm/cc. (a) Vs-Vp trend 
in normalized units. Black line is the fit trend, Eq. ([3|. (b) 
p-Vp trend in normalized units. Trend lines of constant 5, Eq. 
M, axe colored according to the value of fd = fdci.^- Two 
reference points are shown as black dots and labeled. 



C^ = -0.1916, c 
ays — 216 ft/s. 



= 0.017, Ays 

0.421, nc; 



1280 ft/s, Bys = 0.48, 
1.34, n^ ^ 16.4, Hd = 



-5 



ayp = 350 ft/s, A^ = 0.771, B^ 



-3.68 X 10-5 s/ft. 



0.5, fdc = 0.52, Ap = 1.435 gm/cc, Bp = 7.0 x 10 
(s/ft) (gm/cc), Cp = 0.364 gm/cc, and ap = 0.032 gm/cc. 
We have assumed ps = 2.9 gm/cc and pf = 1.0 gm/cc in 
these relationships. We note that there was a four-fold 
decrease in ap by including the Cp term in the regression 
of the datasets. Note the reasonable value of (j>c and the 
critical scaling of i^ = \fjdjjdc- The self similarity of 
the rock structure implied by this model is validated by 
neutron scattering experiments^. Given the value of n^, 
simple arguments indicate that the coordination number 
of the ductile fraction also scales as 0.5. Characteristic 



values of the rock physics parameters are ( = 0.75 ± 0.07 
and £, = 0.70 ± 0.20. Straight forward analysis shows 
that the capture fraction, as defined by Demartini and 
Glinskyl^l, scales as n^/(n<^ — n^) and is approximately 
equal to the reciprocal of the exponent for states away 
from the critical point. This gives a capture fraction of 
92% for this model, and 36% for the floating grain work 
of Demartini and Glinsky-^ . 

We have not explicitly identified the process and there- 
fore the "activation energy" in the definition of ^ = 



1 — exp{—E/EQ). Unlike for the floating grain model, 
changes in the composition are not limited to compaction 
(there is probably a very large amount of diagenesis and 
mineral substitution for shales), and we did not have in- 
formation on what the controlling variables (i.e., effective 
stress or temperature) were for each of the well log sam- 
ples. Practically, this is not a limitation since we are not 
trying to estimate the energy, E, and that value will be 
assumed to be a constant for a stratigraphic layer in our 
analysis. 



III. GEOPHYSICAL FORWARD MODEL 



We now move to developing an understanding of both the P-to-P, Rpp, and the P-to-S, Rps, reflection response for 
an isotropic medium. We start by assuming weak scattering and make both the further assumptions of small contrast 
(that is, Ap/p, Avp/vp, andAvs/vs <C 1) and plane waves (eikonal approximation). The latter is a rather complicated 
assumption on both the frequency and angle of incidence, 0. The expressions for the reflection response will be linear 
in the contrasts due to the flrst approximation, with coefficients that are functions of the angle of incidence, 6, and 
the ratio of the velocities, rs„ = Vg/vp. Expanding to 4th order in 9 leads to the expressions 
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(6) 

(7) 



We have been careful to write these expressions in a bilinear form in terms of the small contrast (i.e., 
Ap/p, Avp/vp,Avs/vs) and the angle of incidence (i.e., 6). This will facilitate the SVD analysis of the next sec- 
tion. Note that the coefficients of this bilinear transformation are only functions of the dimensionless parameter. 



r 



IV. SINGULAR VALUE DECOMPOSITION 

We now move onto understanding the relationships be- 
tween the basic rock physics parameters we wish to know, 
S, and (, and the geophysical measurement. We do this 
by establishing a sequence of linear transformations, then 
examining the important singular value decompositions 
(SVDs) of that compound transformation. The singular 
values will give an understanding of detectability of the 
singular vectors (that is, the required SNR). The singular 
vectors will tell us what views of the measurement to use 
and how they are related to the rock physics. 

We need to understand the decomposition of x^ where 



Mrp, 



R = Me Ma Mrp Ar, 



(9) 



where R = MqA, A^ Ma Ac, and Ac = Mp^p Ar. 

With these definitions in hand we now return to the 
definition of x and write it as 



X^Wd R,n ~WdD Mg Ma Mrp Ar. 



(10) 



We first make an SVD of WaDMg = UiY^iV^ . We de- 
note El as the square matrix formed by dropping the 
zero rows of Si. We also define Ui by dropping the same 
rows of C/i. We make a second SVD of Y^iVi MaMrp = 
C/2S2V2'". It can be shown that 



WdRm-WdDR = x- 



(8) 



R is the theoretical value, Rm is the measured value of 
R, D la a, linear distortion of the measurement of R, 



Wd 



^d'Ud, 



E J is a diagonal matrix of the standard 



deviations, l/cr, and Ud is unitary. Before we continue 
with understanding the linear structure of x^, we need 
to develop the linear structure of R in more detail. It 
can be written as a product of the angle matrix, Mg, 
the reflection matrix, Ma, and the rock physics matrix. 



{U^U,xf{U2U,x)^X^X, 



(11) 



tT^T 



so we are free to redefine x &s U2 U ^ x without changing 
X^ X — X^; that is the probability distribution. There- 
fore, we write 



X 



ul{ulWd)R,, 



Y.2V^Ar. 



(12) 



Let us now make some practical identifications. First of 

— T 

all recognize that U^ Wd transforms Rm into to "stacks" 



where m is the dimension of A. We will denote these 

stacks as Ri. V^ is a 2 x 2 matrix that rotates Ar so 

that they are orthogonal, Ar — V^Ar. Then the m 

J" 

stacks R are projected by t/2 {a. 2 x m matrix) onto the 

two orthogonal rock physics directions. The two singular 
values given by the diagonal matrix E2 give the SNR of 
the estimates of the two orthogonal rock physics param- 
eters. One can directly form the two optimal stacks for 

— T — T 

estimation of the rock physics parameters by t/j U-^ Wd- 
This is still abstract at this point. Let us substitute in 
the rock physics of the shales given in the latter part of 
Sec. [n] Set the multiplicative distortion to the identity 
matrix and the data covariance to a diagonal constant 
of 1. We set the rock physics composition to C — 0.79 
and the geometry to ^ = 0.5. This gives a density of p = 
2.59 gm/cc, compressional velocity of Vp = 14000 ft/s, a 
shear velocity of Vs = 8000 ft/s, a Vp to v^ ratio of Vpg — 
1.75, a Poisson ratio of i^ = 0.26, and a porosity oi (p ~ 
16%. 

For a typical angle of Om ~ 30°, we get the stack 
weights, Ui , shown in Fig. m It should be noted that 
this result is independent of the rock physics, Mrp, and 
the relationship between the rock physics and the A's, 
Ma- In the order of decreasing singular value, or SNR, 
we have Rq the full PP stack, Ri the "full" PS stack (in 
quotes because it is really linearly weighted with 6) , then 
i?2 the AVO PP gradient stack (weighted by 9^ so that it 
is the far offsets minus the near offsets). The series con- 
tinues on with progressively higher order weightings of 
the stacks in an alternating order between the PP and 
the PS data. The singular values scale as Ai ~ 6*^. Con- 
tinuing with the analysis, we analyze the rotation of Ar 
onto an orthogonal system Ar. We find that C is mainly 
the composition variable C and the ^ variable is mainly 

the geometry variable ^. Figurewlshows the U2 transfor- 
mation of the stacks, R, onto the rock physics variables, 
Af. Note that the full PP stack is the main contribu- 
tion to the determination of the composition variable, ^, 
and the "full" PS stack is the main contribution to the 
determination of the geometry variable, ^. The AVO PP 
gradient stack is of minor contribution to either, but it is 
more aligned with ^ and orthogonal to C. The 4th order 
PP, i?4, is totally neghgible. 



V. CONCLUSIONS 

The purpose of this paper has been to establish the un- 
derlying fundamentals of quantitative interpretation for 
unconventional shale reservoirs. This starts with the un- 
derstanding that it is the ductile fraction that controls 
the geomechancial balance between the rocky road joint 
friction of the fractures and the viscous joint friction. 
While this geomechanics is not the subject of this paper, 
otherslJ have found sensitive dependance of the dynamic 
friction to the ductile fraction, and a resulting dramatic 
change in the fracturing efficiency. 




FIG. 2. Stack weights, C/i , as a function of incidence angle, 
0. First set is for PP data, followed by the weights for PS 
data. 
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FIG. 3. Transformation of the stacks onto the rock physics 

— T 

parameters, U2 , for 6m = 30°. 



Inspired by this geomechanical observation, we de- 
veloped and verified, at the mesoscopic level, a rock 
physics model where the three isotropic elastic properties 
are only a function of two parameters, the scaled duc- 
tile fraction, ^ = y/Jd/Jdc, and a composition variable, 
C — 1 — exp {—E/Eq), which captures compaction, diage- 
nesis, and mineral substitution effects. The first variable 
captures changes in the geometric microstructure, that 
is how efficient the rock matrix is in supporting stress - 
modulus per mass or coordination number. The second 
variable captures the compositional properties of the ma- 
trix. It is a remarkable gift of nature that there are only 
two parameters and that one of them is directly related 
to the ductile fraction - the critical parameter for the 
geomechanics. 

The next important question that we answered is how 
this geometry parameter, ^, manifests itself in surface re- 
flection seismic. The equations relating the rock physics 
to the reflectivity were all linearized and an SVD analy- 
sis was done to answer this question. The leading order 
singular value was directly related to the full PP stack 
and the composition, ^. The next order singular value 



was directly related to t he "ful l" PS stack and the ge- 
onietry variable, ^ = yfdiJdc- For reasonable angles 
of reflection, the higher order stacks, which include the 
AVO PP gradient stack, all have SNR which would make 
them hard, if not impossible, to detect. 
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